;;########################################################################
;; anovamoc.lsp
;; Copyright (c) 1993-99 by Forrest W. Young
;; has code for save-model, leverages, ls-means
;;########################################################################



(defmeth anova-model-object-proto :save-model-template (data-object)
"Args: (data-object)
DATA-OBJECT is the object-identification information of a data object. 
The method contains a template for saving the model-object." 
  `(analysis-of-variance    
    :title      ,(send self :title)
    :name       ,(send self :name) 
    :dialog      nil
    :interaction ,(send current-model :interaction)
    :data (data  ,(send data-object :name)
                 :title      ,(send data-object :title)
                 :variables ',(send self :variables)
                 :types     ',(send self :types)
                 :labels    ',(send self :labels)
                 :ways      ',(send self :ways)
                 :classes   ',(send self :classes)
                 :data      ',(send self :data)))
  )

(defmeth anova-model-object-proto :create-data (&key (dialog nil))
  (if (not (eq current-object self)) (setcm self)) 
  (let* ((creator (send *desktop* :selected-icon))
         (data-obj (send self :data-object))
         (nway (send self :nway-model))
         (nways (send self :nways))
         (data-ways (send data-obj :ways))
         (model-ways (send self :ways))
         (response-var (first (send self :variables)))
         (y (send nway :y))
         (x (if data-ways (send nway :x) (send self :class-matrix)))
         (nx (second (size x)))
         (x-var-names (if data-ways (send nway :predictor-names) model-ways))
         (fit-values (send nway :fit-values))
         (raw-residuals (send nway :raw-residuals))
         (residuals (send nway :residuals))
         (studentized-residuals (send nway :studentized-residuals))
         (externally-studentized-residuals 
                               (send nway :externally-studentized-residuals))
         (cooks-distances (send nway :cooks-distances))
         (data (bind-columns y x fit-values residuals raw-residuals 
                    studentized-residuals externally-studentized-residuals 
                    cooks-distances))
         (types (repeat "Numeric" (second (size data))))
         )
    (if (not data-ways)
        (setf (select types (iseq 1 nways)) (repeat "Category" nways)))
    (data (concatenate 'string "Results-" (send self :name))
          :created creator
          :creator-object self
          :title (concatenate 'string "ANOVA Results From " 
                              (send self :title))
          :data (combine data)
          :types types
          :variables (combine response-var x-var-names "FitValues" 
                              "Residuals" "RawResids"  
                              "StudResids" "ExtStudResids" 
                              "CooksDists")
          :labels 
          (if (not data-ways)
              (send self :labels)
              (repeat (send data-obj :labels) (send data-obj :cellfreqs)))
          )))


(defmeth anova-model-object-proto :lev (sel)
  (let* ((nway (send self :nway-model))
         (ncols (select (array-dimensions (send nway :x)) 1))
         (ident-mat (identity-matrix ncols))
         (nrows nil)
         (mats (mapcar #'array-dimensions (send nway :design-matrices)))
         (x nil)
         (g nil)
         (i 0)
         (j 0)
         (levmat nil)
         (mat nil)
         )
    (dotimes (i (length mats))
             (setf x (remove nil (combine x 
                     (list (select (select mats i) 1))))))
    (when (not (= sel 0))
          (setf nrows (select x (1- sel)))
          (setf g (if (not (= sel 1))
                      (apply #'+ (select x (iseq (1- sel))))))
          (if (= sel 1) (setf g 0))
          (setf levmat (select ident-mat (iseq g (+ g (1- nrows))) (iseq ncols)))
          (setf mat (send self :sall-lev levmat)))
    (when (= sel 0)
          (setf mat (list (send nway :fit-values) (send nway :y))))
    mat))

(defmeth anova-model-object-proto :sall-lev (L)
"Args: L
L is an indicator matrix for a source in the anova design. L has a column for every term in the design. L has one row for each level of a source, a row has a single 1 indicating the factors position in the design."
  (let* ((nway    (send self :nway-model))
         (x       (send nway :x))
         (xpxiLp  (matmult (inverse (matmult (transpose x) x))(transpose L)))
         (xxpxiLp (matmult x xpxiLp))
         (b       (bind-columns (rest (send nway :coef-estimates))))
         (lambda  (matmult (inverse (matmult L xpxiLp)) L b))
         (vx      (+ (mean (send nway :y)) 
                     (combine (matmult xxpxiLp lambda))))
         (vy      (+ (send nway :residuals) vx))
         )
    (list vx vy)))

(defmeth anova-model-object-proto :make-leverages ()
"Args: none
Returns an nobs x nsources matrix of Sall's Leverage values for the predictor source variables."
  (let* ((levs nil)
         (nway-model (send self :nway-model))
         (n-sources 
          (- (length (send nway-model :source-degrees-of-freedom)) 2))
         (n-obs (send nway-model :nobs))
         )
    (dotimes (i n-sources)
             (setf levs (combine levs (first (send self :lev (1+ i))))))
    (transpose (matrix (list n-sources n-obs) (rest levs)))))

(defmeth anova-model-object-proto :grouped-values (values model sel)
  (let* ((indicator 
          (nth sel (send self :indicator-matrices)))
         (nclasses (second (size indicator)))
         (members nil)
         (grouped-values nil))
    (dotimes (i nclasses)
             (setf members (select values (which (= 1 (col indicator i)))))
             (setf grouped-values (append grouped-values (list members))))
    grouped-values))

(defmeth anova-model-object-proto :make-lsmeans ()
  (let* ((grouped-data)
         (ls-means)
         (point-labels (send self :labels))
         (num-sources (send self :nsources)))
    (dotimes (i num-sources)
             (setf grouped-data 
                   (first (send self :grouped-data i point-labels)))
             (setf ls-means (append ls-means (list (mapcar #'mean grouped-data))))
             )
    ls-means))
